package ec.app.regression;
import ec.util.*;
import ec.*;
import ec.gp.*;
import ec.gp.koza.*;
import ec.simple.*;
import java.io.*;
import java.util.*;

/* 
 * Benchmarks.java
 * U    
 * Created: Thu Jul 14 10:35:11 EDT 2011
 * By: Sean Luke
 *
 * This is an expansion of the Regression.java file to provide a first cut at a standardized
 * regression benchmark for Genetic Programming.  The package provides both training and
 * testing data and functions for the following problems:
 
 <ol>
 <li>John Koza's three  problems (quartic, quintic, and sextic, drawn from his GP-1 and GP-2 books.  These are known as <b>koza-1</b> through <b>koza-3</b>
 <li>Twelve problems drawn from "Semantically-based Crossover in Genetic Programming: Application to Real-valued Symbolic Regression" (GPEM 2011),
 by Nguyen Quang Uy, Nguyen Xuan Hoai, Michael O’Neill, R.I. McKay, Edgar Galv ́an-L ́opez.  These are known as
 <b>nguyen-1</b> through <b>nguyen-10</b>
 <li>Fifteen problems drawn from "Accuracy in Symbolic Regression" (GPEM 2011), by Michael Korns.  These are known as <b>KORNS1</b> through <b>KORNS15</b>
 <li>Fifteen problems drawn from "Improving Symbolic Regression with Interval Arithmetic and Linear Scaling" (EuroGP 2003), by Maarten Keijzer.  These are known as <b>keijzer-1</b> through <b>keijzer-15</b>
 <li>Fifteen problems drawn from "Order of Nonlinearity as a Complexity Measure for Models Generated by Symbolic Regression via Pareto Genetic Programming" (IEEE TransEC 13:2), by Ekaterina J. Vladislavleva, Guido F. Smits, and Dick den Hertog.  These are known as <b>vladislavleva-1</b> through <b>vladislavleva-8</b>
 <li>Two problems drawn from "Evolutionary Consequences of Coevolving Targets" (Evolutionary Computation 5(4)) by Ludo Pagie and Paulien Hogeweg, 1997.
 <li>You can also provide your own data sets via a text file.

 <p>These problems differ in a variety of ways.  Some have both training and testing data points.  Some have one variable, others up to five variables.  Some build their random points based on random samples, other based on a grid of samples distributed through the space.  Different problems also use different function sets.
 
 <p>The functions below are all described precisely in the paper "Genetic Programming Needs Better Benchmarks" (GECCO 2012) by James McDermmott, David R. white, Sean Luke, Luca Manzoni, Mauro Castelli, Leonardo Vanneschi, Wojciech Jaskowski, Krzysztof Krawiec, Robin Harper, Kenneth De Jong, and Una-May O'Reilly.  The descriptions, shown in Tables 2 and 3, explain the function set, number of variables, objective function, training set, testing set, and source of each problem.
 
 <p>In addition we include one more function: PAGIE2, a 3-variable version of PAGIE1. We describe PAGIE2 as follows:
 
 <ul>
 <li>Variables: 3
 li>Function: (1 / (1 + x^(-4)) + 1 / (1 + y^(-4)) + 1 / (1 + z^(-4)))
 <li>Testing Set: Grid points from -5 to 5 respectively in each dimension, spaced out in intervals of 0.4.
 <li>Training Set: none
 <li>Function Set: standard Koza with three terminals (X1, X2, X3).
 </ul>

 *
 *
 */

/**
 * Benchmarks by various people in the literature.
 *
 */
public class Benchmarks extends GPProblem implements SimpleProblemForm
    {
    private static final long serialVersionUID = 1;

    //////////////////////////// problem tags
        
    // From GP-1 and GP-2
    public static final int KOZA1 = 0 ;
    public static final int KOZA2 = 1 ;
    public static final int KOZA3 = 2 ;
        
    // From Semantically-based Crossover in Genetic Programming: Application to Real-valued Symbolic Regression
    // Nguyen Quang Uy, Nguyen Xuan Hoai, Michael O’Neill, R.I. McKay, Edgar Galv ́an-L ́opez
    // GPEM 2011
    public static final int NGUYEN1 = 3 ;
    public static final int NGUYEN2 = 4 ;  // identical to standard Koza
    public static final int NGUYEN3 = 5 ;
    public static final int NGUYEN4 = 6 ;
    public static final int NGUYEN5 = 7 ;
    public static final int NGUYEN6 = 8 ;
    public static final int NGUYEN7 = 9 ;
    public static final int NGUYEN8 = 10 ;
    public static final int NGUYEN9 = 11 ;
    public static final int NGUYEN10 = 12 ;
    
    // From Pagie and Hogeweg 1997: Ludo Pagie, Paulien Hogeweg: Evolutionary Consequences of Coevolving
    // Targets. Evolutionary Computation 5(4): 401-418 (1997).
    public static final int PAGIE1 = 13;        // 2D
    public static final int PAGIE2 = 14;        // 3D
        
    // From Accuracy in Symbolic Regression -- GPTP 2011
    public static final int KORNS1 = 15 ;
    public static final int KORNS2 = 16 ;
    public static final int KORNS3 = 17 ;
    public static final int KORNS4 = 18 ;
    public static final int KORNS5 = 19 ;
    public static final int KORNS6 = 20 ;
    public static final int KORNS7 = 21 ;
    public static final int KORNS8 = 22 ;
    public static final int KORNS9 = 23 ;
    public static final int KORNS10 = 24 ;
    public static final int KORNS11 = 25 ;
    public static final int KORNS12 = 26 ;
    public static final int KORNS13 = 27 ;
    public static final int KORNS14 = 28 ;
    public static final int KORNS15 = 29 ;
        
    // From Improving Symbolic Regression with Interval Arithmetic and Linear Scaling
    // Maarten Keijzer -- EuroGP 2003
    public static final int KEIJZER1 = 30 ;
    public static final int KEIJZER2 = 31 ;
    public static final int KEIJZER3 = 32 ;
    public static final int KEIJZER4 = 33 ;
    public static final int KEIJZER5 = 34 ;
    public static final int KEIJZER6 = 35 ;
    public static final int KEIJZER7 = 36 ;                         
    public static final int KEIJZER8 = 37 ;
    public static final int KEIJZER9 = 38 ;
    public static final int KEIJZER10 = 39 ;
    public static final int KEIJZER11 = 40 ;
    public static final int KEIJZER12 = 41 ;
    public static final int KEIJZER13 = 42 ;
    public static final int KEIJZER14 = 43 ;
    public static final int KEIJZER15 = 44 ;
        
    // From Order of Nonlinearity as a Complexity Measure for Models Generated by Symbolic Regression via Pareto Genetic Programming
    // Ekaterina J. Vladislavleva, Guido F. Smits, and Dick den Hertog
    public static final int VLADISLAVLEVA1 = 45 ;
    public static final int VLADISLAVLEVA2 = 46 ;
    public static final int VLADISLAVLEVA3 = 47 ;
    public static final int VLADISLAVLEVA4 = 48 ;
    public static final int VLADISLAVLEVA5 = 49 ;
    public static final int VLADISLAVLEVA6 = 50 ;                           
    public static final int VLADISLAVLEVA7 = 51 ;
    public static final int VLADISLAVLEVA8 = 52 ;

        
        
        
    // parameter names
    public static final String names[] = 
        { 
        "koza-1", "koza-2", "koza-3",
        "nguyen-1", "nguyen-2", "nguyen-3", "nguyen-4", "nguyen-5", "nguyen-6", "nguyen-7", "nguyen-8", "nguyen-9", "nguyen-10",
        "pagie-1", "pagie-2",
        "korns-1", "korns-2", "korns-3", "korns-4", "korns-5", "korns-6", "korns-7", "korns-8", "korns-9", "korns-10", "korns-11", "korns-12", "korns-13", "korns-14", "korns-15",
        "keijzer-1", "keijzer-2", "keijzer-3", "keijzer-4", "keijzer-5", "keijzer-6", "keijzer-7", "keijzer-8", "keijzer-9", "keijzer-10", "keijzer-11", "keijzer-12", "keijzer-13", "keijzer-14", "keijzer-15", 
        "vladislavleva-1", "vladislavleva-2", "vladislavleva-3", "vladislavleva-4", "vladislavleva-5", "vladislavleva-6", "vladislavleva-7", "vladislavleva-8"
        }; 
                
                
    // expected function sets.  "fn" means "function set with n terminals x_1 ... x_n"
    public static final String fs[] = 
        {
        "koza1", "koza1", "koza1",
        "koza1", "koza1", "koza1", "koza1", "koza1", "koza1", "koza1", "koza1", "koza2", "koza2",
        "koza2", "koza3",
        "korns5", "korns5", "korns5", "korns5", "korns5", "korns5", "korns5", "korns5", "korns5", "korns5", "korns5", "korns5", "korns5", "korns5", "korns5", 
        "keijzer1", "keijzer1", "keijzer1"/* revised from GECCO paper */,  "keijzer1", "keijzer3", "keijzer1", "keijzer1", "keijzer1", "keijzer1", "keijzer2", "keijzer2", "keijzer2", "keijzer2", "keijzer2", "keijzer2",
        "vladislavleva-b2", "vladislavleva-c1", "vladislavleva-c2", "vladislavleva-a5", "vladislavleva-a3", "vladislavleva-b2", "vladislavleva-c2", "vladislavleva-a2"
        };
        
        
    
    // function sets with various variable lengths
    public static final String fs_vars[][] = 
        {
        { },
        { "koza1", "keijzer1", "vladislavleva-c1" },
        { "koza2", "keijzer2", "vladislavleva-a2", "vladislavleva-b2", "vladislavleva-c2" },
        { "koza3", "keijzer3", "vladislavleva-a3" },
        { },
        { "korns5", "vladislavleva-a5" },
        };
        
        
    /** Hyperbolic Arc Sin -- not standard in Java Math library */
    static double asinh(double x)
        {
        return Math.log(x + Math.sqrt(x*x + 1.0));
        }
        
        
    //// Input Sample Generation
        
        
    /** Produce random sample points between min and max, inclusive, in each dimension.  */
    public double[][] generateRandomSamples(EvolutionState state, double[] min, double[] max, int numPoints, int threadnum)
        {
        int vars = max.length;
        double[][] d = new double[numPoints][vars];
        for(int i = 0 ; i < d.length; i++)
            for(int j = 0; j < vars; j++)
                // below we're doing nextDouble(true, true), which means to select from the
                // FULLY CLOSED interval [0.0, 1.0], including both 0.0 and 1.0.
                d[i][j] = state.random[threadnum].nextDouble(true, true) * (max[j] - min[j]) + min[j];
        return d;
        }
                
    /** Produce random sample points between min and max, inclusive, in one dimension.  */
    public double[][] generateRandomSamples(EvolutionState state, double min, double max, int numPoints, int threadnum)
        {
        return generateRandomSamples(state, new double[] { min }, new double[] { max }, numPoints, threadnum);
        }

    // recursive trick to dump the full mesh into a bag.  Enter this by setting variable to 0,  Yuck, expensive.  But O(n).
    void buildIntervalPoints(EvolutionState state, ArrayList list, double[] min, double[] max, double[] interval, double current[], int variable, int threadnum)
        {
        if (variable == min.length)  // we're out of variables, base case
            {
            double[] d = new double[min.length];
            for(int i = 0; i < d.length; i++)
                d[i] = current[i];      // not sure if System.arraycopy would be faster, probably not in this case
            list.add(d);
            }
        else
            {
            int jumps = (int)((max[variable] - min[variable]) / interval[variable]) + 1;

            for(int j = 0; j < jumps; j++)
                {
                current[variable] = min[variable] + interval[variable] * j;
                buildIntervalPoints(state, list, min, max, interval, current, variable + 1, threadnum);
                }
            }
        }



    /** Produce sample points evenly spaced out between min and max in each dimension, with the given spacing interval per-dimension.  */
    public double[][] generateIntervalSpacedSamples(EvolutionState state, double[] min, double[] max, double[] interval, int threadnum)
        {
        // gather all the points in the mesh recursively
        ArrayList list = new ArrayList();
        double[] current = new double[min.length];
        buildIntervalPoints(state, list, min, max, interval, current, 0, threadnum);
                
        // Convert to an array
        return (double[][])(list.toArray(new double[0][]));
        }
                
                
    /** Produce sample points evenly spaced out between min and max in one dimension, with the given spacing interval.  One dimension only. */
    public double[][] generateIntervalSpacedSamples(EvolutionState state, double min, double max, double interval, int threadnum)
        {
        return generateIntervalSpacedSamples(state, new double[] { min }, new double[] { max }, new double[] { interval }, threadnum);
        }
                
        
    /** Produce sample points for a given benchmark problem.  */
    public double[][] trainPoints(EvolutionState state, int benchmark, int threadnum)
        {
        switch(benchmark)
            {
            case KOZA1:
            case KOZA2:
            case KOZA3:
                return generateRandomSamples(state, -1, 1, 20, threadnum);
                                
            case NGUYEN1:
            case NGUYEN2:
            case NGUYEN3:
            case NGUYEN4:
            case NGUYEN5:
            case NGUYEN6:
                return generateRandomSamples(state, -1, 1, 20, threadnum);
                        
            case NGUYEN7:
                return generateRandomSamples(state, 0, 2, 20, threadnum);

            case NGUYEN8:
                return generateRandomSamples(state, 0, 4, 20, threadnum);

            case NGUYEN9:
            case NGUYEN10:
                return generateRandomSamples(state, new double[] { 0, 0 }, new double[] { 1, 1 }, 100, threadnum);
                        
            case PAGIE1:
                return generateIntervalSpacedSamples(state, new double[] { -5.0, -5.0 }, new double[] { 5.0, 5.0 }, new double[] { 0.4, 0.4 }, threadnum);
            case PAGIE2:
                return generateIntervalSpacedSamples(state, new double[] { -5.0, -5.0, -5.0 }, new double[] { 5.0, 5.0, 5.0 }, new double[] { 0.4, 0.4, 0.4 }, threadnum);
            
            
            case KORNS1:
            case KORNS2:
            case KORNS3:
            case KORNS4:
            case KORNS5:
            case KORNS6:
            case KORNS7:
            case KORNS8:
            case KORNS9:
            case KORNS10:
            case KORNS11:
            case KORNS12:
            case KORNS13:
            case KORNS14:
            case KORNS15:
                return generateRandomSamples(state, new double[] { -50, -50, -50, -50, -50 }, new double[] { 50, 50, 50, 50, 50 }, 10000, threadnum);  // 10000 !!
                        
            case KEIJZER1:
                return generateIntervalSpacedSamples(state, -1, 1, 0.1, threadnum);
                                
            case KEIJZER2:
                return generateIntervalSpacedSamples(state, -2, 2, 0.1, threadnum);

            case KEIJZER3:
                return generateIntervalSpacedSamples(state, -3, 3, 0.1, threadnum);

            case KEIJZER4:
                return generateIntervalSpacedSamples(state, 0, 10, 0.05, threadnum);
                                
            case KEIJZER5:      // revised from GECCO paper
                return generateRandomSamples(state, new double[] { -1, 1, -1 }, new double[] { 1, 2, 1 }, 1000, threadnum);
                                
            case KEIJZER6:
                return generateIntervalSpacedSamples(state, 1, 50, 1, threadnum);
                                
            case KEIJZER7:                          
                return generateIntervalSpacedSamples(state, 1, 100, 1, threadnum);
            case KEIJZER8:
            case KEIJZER9:
                return generateIntervalSpacedSamples(state, 0, 100, 1, threadnum);

            case KEIJZER10:
                return generateRandomSamples(state, new double[] { 0, 0 }, new double[] { 1, 1 }, 100, threadnum);

            case KEIJZER11:
            case KEIJZER12:
            case KEIJZER13:
            case KEIJZER14:
            case KEIJZER15:
                return generateRandomSamples(state, new double[] { -3, -3 }, new double[] { 3, 3 }, 20, threadnum);
                                
            case VLADISLAVLEVA1:
                return generateRandomSamples(state, new double[] { 0.3, 0.3 }, new double[] { 4, 4 }, 100, threadnum);
                                
            case VLADISLAVLEVA2:      // revised from GECCO paper
                return generateIntervalSpacedSamples(state, 0.05, 10, 0.1, threadnum);
                                
            case VLADISLAVLEVA3:
                return generateIntervalSpacedSamples(state, new double[] { 0.05, 0.05 }, new double[] { 10, 10.05 }, new double[] { 0.1, 2 }, threadnum);
                                
            case VLADISLAVLEVA4:
                return generateRandomSamples(state, new double[] { 0.05, 0.05, 0.05, 0.05, 0.05}, new double[] { 6.05, 6.05, 6.05, 6.05, 6.05 }, 1024, threadnum);
                                
            case VLADISLAVLEVA5:
                return generateRandomSamples(state, new double[] { 0.05, 1, 0.05 }, new double[] { 2, 2, 2 }, 300, threadnum);
                                
            case VLADISLAVLEVA6:                            
                return generateRandomSamples(state, new double[] { 0.1, 0.1 }, new double[] { 5.9, 5.9 }, 30, threadnum);

            case VLADISLAVLEVA7:
                return generateRandomSamples(state, new double[] { 0.05, 0.05 }, new double[] { 6.05, 6.05 }, 300, threadnum);

            case VLADISLAVLEVA8:
                return generateRandomSamples(state, new double[] { 0.05, 0.05 }, new double[] { 6.05, 6.05 }, 50, threadnum);
                        
            default:
                return null;
            }
        }



    /** Produce test sample points for a given benchmark problem, to test generalization.  */
    public double[][] testPoints(EvolutionState state, int benchmark, int threadnum, double[][] trainpoints)
        {
        switch(benchmark)
            {
            case KOZA1:
            case KOZA2:
            case KOZA3:
                return trainpoints;
                                
            case NGUYEN1:
            case NGUYEN2:
            case NGUYEN3:
            case NGUYEN4:
            case NGUYEN5:
            case NGUYEN6:
            case NGUYEN7:
            case NGUYEN8:
            case NGUYEN9:
            case NGUYEN10:
                return trainpoints;
                
            case PAGIE1:
            case PAGIE2:
                return trainpoints;
                        
            case KORNS1:
            case KORNS2:
            case KORNS3:
            case KORNS4:
            case KORNS5:
            case KORNS6:
            case KORNS7:
            case KORNS8:
            case KORNS9:
            case KORNS10:
            case KORNS11:
            case KORNS12:
            case KORNS13:
            case KORNS14:
            case KORNS15:
                return generateRandomSamples(state, new double[] { -50, -50, -50, -50, -50 }, new double[] { 50, 50, 50, 50, 50 }, 10000, threadnum);  // 10000 !!
                        
            case KEIJZER1:
                return generateIntervalSpacedSamples(state, -1, 1, 0.001, threadnum);
                                
            case KEIJZER2:
                return generateIntervalSpacedSamples(state, -2, 2, 0.001, threadnum);

            case KEIJZER3:
                return generateIntervalSpacedSamples(state, -3, 3, 0.001, threadnum);

            case KEIJZER4:
                return generateIntervalSpacedSamples(state, 0.05, 10.05, 0.05, threadnum);
                                
            case KEIJZER5:      // revised from GECCO paper
                return generateRandomSamples(state, new double[] { -1, 1, -1 }, new double[] { 1, 2, 1 }, 10000, threadnum);  // 10000 cases for testing, different than for training
                                
            case KEIJZER6:
                return generateIntervalSpacedSamples(state, 1, 120, 1, threadnum);
                                
            case KEIJZER7:                          
                return generateIntervalSpacedSamples(state, 1, 100, 0.1, threadnum);
            case KEIJZER8:
            case KEIJZER9:
                return generateIntervalSpacedSamples(state, 0, 100, 0.1, threadnum);

            case KEIJZER10:
                return generateIntervalSpacedSamples(state, new double[] { 0, 0 }, new double[] { 1, 1 }, new double[] { 0.01, 0.01 }, threadnum);

            case KEIJZER11:
            case KEIJZER12:
            case KEIJZER13:
            case KEIJZER14:
            case KEIJZER15:
                return generateIntervalSpacedSamples(state, new double[] { -3.0, -3.0 }, new double[] { 3.0, 3.0 }, new double[] { 0.01, 0.01 }, threadnum);
                                
            case VLADISLAVLEVA1:      // revised from GECCO paper
                return generateIntervalSpacedSamples(state, new double[] { -0.2, -0.2 }, new double[] { 4.2, 4.2 }, new double[] { 0.1, 0.1 }, threadnum);
                                
            case VLADISLAVLEVA2:
                return generateIntervalSpacedSamples(state, -0.5, 10.5, 0.05, threadnum);
                                
            case VLADISLAVLEVA3:
                return generateIntervalSpacedSamples(state, new double[] { -0.5, -0.5 }, new double[] { 10.5, 10.5 }, new double[] { 0.05, 0.5 }, threadnum);  // note 0.05 and 0.5, and also 10.5 which is different from training
                                
            case VLADISLAVLEVA4:
                return generateRandomSamples(state, new double[] { -0.25, -0.25, -0.25, -0.25, -0.25}, new double[] { 6.35, 6.35, 6.35, 6.35, 6.35 }, 5000, threadnum);
                                
            case VLADISLAVLEVA5:
                return generateIntervalSpacedSamples(state, new double[] { -0.05,  0.95, -0.05 }, new double[] { 2.1, 2.05, 2.1 }, new double[] { 0.15, 0.15, 0.1 }, threadnum);  // note 0.05 and 0.5, and also 10.5 which is different from training
                                
            case VLADISLAVLEVA6:                            
                return generateIntervalSpacedSamples(state, new double[] { -0.05, -0.05 }, new double[] { 6.05, 6.05 }, new double[] { 0.02, 0.02 }, threadnum);  // note 0.05 and 0.5, and also 10.5 which is different from training

            case VLADISLAVLEVA7:
                return generateRandomSamples(state, new double[] { -0.25, -0.25 }, new double[] { 6.35, 6.35 }, 1000, threadnum);

            case VLADISLAVLEVA8:
                return generateIntervalSpacedSamples(state, new double[] { -0.25, -0.25 }, new double[] { 6.35, 6.35 }, new double[] { 0.2, 0.2 }, threadnum);  // note 0.05 and 0.5, and also 10.5 which is different from training

            default:
                return null;
            }
        }







    /** Return the function applied to the given data by benchmark problem.  */
    public double func(EvolutionState state, double[] xs, int benchmark) throws IllegalArgumentException
        {
        double x = xs[0];
        double y = (xs.length > 1 ? xs[1] : 0);
        double z = (xs.length > 2 ? xs[2] : 0);
                
        switch(benchmark)
            {
            case KOZA1:
                return x*x*x*x + x*x*x + x*x + x;       // traditional
            case KOZA2:
                return x*x*x*x*x - 2.0*x*x*x + x;       // Quintic, from  J. R. Koza, GP II, 1994
            case KOZA3:
                return x*x*x*x*x*x - 2.0*x*x*x*x + x*x;  // Sextic, from J. R. Koza, GP II, 1994
            case NGUYEN1:
                return x*x*x + x*x + x;
            case NGUYEN2:                                                   // identical to KOZA1
                return x*x*x*x + x*x*x + x*x + x;
            case NGUYEN3:
                return x*x*x*x*x + x*x*x*x + x*x*x + x*x + x;
            case NGUYEN4:
                return x*x*x*x*x*x + x*x*x*x*x + x*x*x*x + x*x*x + x*x + x;
            case NGUYEN5:
                return Math.sin(x*x) * Math.cos(x) - 1.0;
            case NGUYEN6:
                return Math.sin(x) + Math.sin(x*x + x);
            case NGUYEN7:
                return Math.log(x+1) + Math.log(x*x + 1.0);
            case NGUYEN8:                           // Note this presumes you don't have sqrt(x) in your function set!
                return Math.sqrt(x);
            case NGUYEN9:
                return Math.sin(x) + Math.sin(y*y);
            case NGUYEN10:
                return 2 * Math.sin(x) * Math.cos(y);
                //case NGUYEN11:
                //    return Math.pow(x, y);
                //case NGUYEN12:
                //    return x*x*x*x - x*x*x + (y*y)/2.0 - y;
            case PAGIE1:
                // otherwise known as 1 / (1 + Math.pow(x, -4)) + 1 / (1 + Math.pow(y, -4))
                return 1.0 / (1.0 + 1.0 / (x * x * x * x)) + 1.0 / (1.0 + 1.0 / (y * y * y * y)); 
            case PAGIE2:
                // otherwise known as (1 / (1 + Math.pow(x, -4)) + 1 / (1 + Math.pow(y, -4)) + 1 / (1 + Math.pow(z, -4)));
                return 1.0 / (1.0 + 1.0 / ( x * x * x * x)) + 1.0 / (1.0 + 1.0 / (y * y * y * y)) + 1.0 / (1.0 + 1.0 / (z * z * z * z));
            case KORNS1:
                return 1.57 + (24.3 * xs[3]);
            case KORNS2:
                return 0.23 + (14.2 * ((xs[3]+ xs[1])/(3.0 * xs[4])));
            case KORNS3:
                return -5.41 + (4.9 * (((xs[3] - xs[0]) + (xs[1]/xs[4])) / (3 * xs[4])));
            case KORNS4:
                return -2.3 + (0.13 * Math.sin(xs[2]));
            case KORNS5:
                return 3.0 + (2.13 * Math.log(xs[4]));
            case KORNS6:
                return 1.3 + (0.13 * Math.sqrt(xs[0]));
            case KORNS7:
                return 213.80940889 - (213.80940889 * Math.exp(-0.54723748542 * xs[0]));
            case KORNS8:
                return 6.87 + (11.0 * Math.sqrt(7.23 * xs[0] * xs[3] * xs[4]));
            case KORNS9:        // revised from GECCO paper
                return Math.sqrt(xs[0]) / Math.log(xs[1]) * Math.exp(xs[2]) / (xs[3] * xs[3]);
            case KORNS10:
                return 0.81 + (24.3 * (((2.0 * xs[1]) + (3.0 * (xs[2] * xs[2]))) / ((4.0 * (xs[3]*xs[3]*xs[3])) + (5.0 * (xs[4]*xs[4]*xs[4]*xs[4])))));
            case KORNS11:
                return 6.87 + (11.0 * Math.cos(7.23 * xs[0]*xs[0]*xs[0]));
            case KORNS12:
                return 2.0 - (2.1 * (Math.cos(9.8 * xs[0]) * Math.sin(1.3 * xs[4])));
            case KORNS13:
                return 32.0 - (3.0 * ((Math.tan(xs[0]) / Math.tan(xs[1])) * (Math.tan(xs[2])/Math.tan(xs[3]))));
            case KORNS14:
                return 22.0 - (4.2 * ((Math.cos(xs[0]) - Math.tan(xs[1]))*(Math.tanh(xs[2])/Math.sin(xs[3]))));
            case KORNS15:
                return 12.0 - (6.0 * ((Math.tan(xs[0])/Math.exp(xs[1])) * (Math.log(xs[2]) - Math.tan(xs[3]))));
            case KEIJZER1:  // fall thru
            case KEIJZER2:  // fall thru
            case KEIJZER3:
                return 0.3 * x * Math.sin(2 * Math.PI * x);
            case KEIJZER4:
                return x*x*x * Math.exp(-x)*Math.cos(x)*Math.sin(x)* (Math.sin(x)*Math.sin(x)*Math.cos(x) - 1);
            case KEIJZER5:
                return (30.0 * x * z) / ((x - 10.0) * y * y);
            case KEIJZER6:
                { 
                double sum = 0;
                double fx = Math.floor(x);
                for(int i = 1; i < fx + 1; i++)  // up to and including floor(x)
                    sum += (1.0 / i);
                return sum; 
                }
            case KEIJZER7:                          // Note this presumes you don't have log(x) in your function set!
                return Math.log(x);
            case KEIJZER8:                          // same as NGUYEN8
                return Math.sqrt(x);
            case KEIJZER9:
                return asinh(x);   // Not a function in Math
            case KEIJZER10:
                return Math.pow(x, y);
            case KEIJZER11:
                return x * y + Math.sin((x - 1.0) * (y - 1.0));
            case KEIJZER12:
                return x*x*x*x - x*x*x + y*y/2.0 - y;
            case KEIJZER13:
                return 6.0 * Math.sin(x) * Math.cos(y);
            case KEIJZER14:
                return 8.0 / (2.0 + x*x + y*y);
            case KEIJZER15:
                return x*x*x / 5.0 + y*y*y/2.0 - y - x;
            case VLADISLAVLEVA1:
                return Math.exp(-(x-1)*(x-1)) / (1.2 + (y - 2.5)*(y-2.5));
            case VLADISLAVLEVA2:
                return Math.exp(-x)*x*x*x*Math.cos(x)*Math.sin(x)*(Math.cos(x)*Math.sin(x)*Math.sin(x) - 1);
            case VLADISLAVLEVA3:
                return Math.exp(-x)*x*x*x*Math.cos(x)*Math.sin(x)*(Math.cos(x)*Math.sin(x)*Math.sin(x) - 1) * (y - 5);
            case VLADISLAVLEVA4:
                {
                double sum = 0;
                for(int i = 0; i < 5; i++)
                    sum += (xs[i] - 3) * (xs[i] - 3);
                return 10.0 / (5.0 + sum);
                }
            case VLADISLAVLEVA5:
                return (30.0 * (x - 1.0) * (z - 1.0)) / (y * y * (x - 10.0));
            case VLADISLAVLEVA6:
                return 6.0 * Math.sin(x) * Math.cos(y);
            case VLADISLAVLEVA7:
                return (x - 3.0) * (y - 3.0) + 2 * Math.sin((x - 4.0) * (y - 4.0));
            case VLADISLAVLEVA8:
                return ((x - 3.0) * (x - 3.0) * (x - 3.0) * (x - 3.0) + (y - 3.0) * (y - 3.0) * (y - 3.0) - (y - 3.0)) / ((y - 2.0) * (y - 2.0) * (y - 2.0) * (y - 2.0) + 10.0);
            default:
                throw new IllegalArgumentException("Invalid benchmark value " + benchmark);
            }
        // never reaches here
        }
        




    /////// Computation of delta error between expected value and provided value



    final static double PROBABLY_ZERO = 1.11E-15;
    final static double BIG_NUMBER = 1.0e15;                // the same as lilgp uses

    /** Returns the error between the result and the expected result of a single
        data point. */
    public double error(double result, double expectedResult)
        {
        double delta = Math.abs(result - expectedResult);

        // It's possible to get NaN because cos(infinity) and
        // sin(infinity) are undefined (hence cos(exp(3000)) zings ya!)
        // So since NaN is NOT =,<,>,etc. any other number, including
        // NaN, we're CAREFULLY wording our cutoff to include NaN.

        if (! (delta < BIG_NUMBER ) )   // *NOT* (delta >= BIG_NUMBER)
            delta = BIG_NUMBER;

        // very slight math errors can creep in when evaluating
        // two equivalent by differently-ordered functions, like
        // x * (x*x*x + x*x)  vs. x*x*x*x + x*x
        // So we're assuming that very small values are actually zero

        else if (delta < PROBABLY_ZERO)  // slightly off
            delta = 0.0;
        return delta;
        }









    ///// Setup

    // parameters
    public static final String P_TESTING_FILE = "testing-file";
    public static final String P_TRAINING_FILE = "training-file";
    public static final String P_PROBLEM_TYPE = "type";

    public double[] currentValue;
    
    // these are read-only during evaluation-time, so
    // they can be just light-cloned and not deep cloned.
    // cool, huh?
    
    public double[][] trainingInputs;
    public double[] trainingOutputs;
    public double[][] testingInputs;
    public double[] testingOutputs;

    // don't bother cloning the inputs and outputs; they're read-only :-)
    // don't bother cloning the current value, it's only set during evaluation

    public void setup(EvolutionState state, Parameter base)
        {
        // very important, remember this
        super.setup(state,base);

        // verify our input is the right class (or subclasses from it)
        if (!(input instanceof RegressionData))
            state.output.fatal("GPData class must subclass from " + RegressionData.class,
                base.push(P_DATA), null);

        // should we load our x parameters from a file, or generate them randomly?
        InputStream training_file = state.parameters.getResource(base.push(P_TRAINING_FILE), null);
        InputStream testing_file = state.parameters.getResource(base.push(P_TESTING_FILE), null);
        String problem = state.parameters.getString(base.push(P_PROBLEM_TYPE), null);
        int benchmark = -1;
                
        if (problem == null)
            {
            state.output.message("Loading benchmark data from files");
            if ((testing_file == null || training_file == null))            // must provide both
                {
                state.output.fatal("If you don't specify a problem type, you must provide a training file and a testing file",
                    (training_file == null ? base.push(P_TRAINING_FILE) : base.push(P_TESTING_FILE)));
                }
            else  // load from files
                {
                try
                    {
                    int numInputs = 0;
                                        
                    // first load the number of input variables
                    Scanner scan = new Scanner(training_file);
                    if (scan.hasNextInt()) 
                        numInputs = scan.nextInt();
                    else state.output.fatal("Number of input variables not provided at beginning of training file ", base.push(P_TRAINING_FILE), null);
                                        
                    // Load into an array list each element
                    ArrayList input = new ArrayList();
                    ArrayList output = new ArrayList();
                    while(scan.hasNextDouble())
                        {
                        double[] in = new double[numInputs];
                        double out = 0;
                        for(int i = 0; i < numInputs; i++)
                            {
                            if (scan.hasNextDouble())
                                in[i] = scan.nextDouble();
                            else state.output.fatal("Non-normal number of data points in training file ", base.push(P_TRAINING_FILE), null);
                            }
                        if (scan.hasNextDouble())
                            out = scan.nextDouble();
                        else state.output.fatal("Non-normal number of data points in training file ", base.push(P_TRAINING_FILE), null);
                        input.add(in);
                        output.add(new Double(out));
                        }
                                        
                    // dump to arrays
                    int len = input.size();
                    trainingInputs = new double[len][numInputs];
                    trainingOutputs = new double[len];
                    for(int i = 0; i < len; i++)
                        {
                        trainingInputs[i] = (double[])(input.get(i));
                        trainingOutputs[i] = ((Double)(output.get(i))).doubleValue();
                        }

        
                    // same thing for testing


                    scan = new Scanner(testing_file);
                    if (scan.hasNextInt()) 
                        numInputs = scan.nextInt();
                    else state.output.fatal("Number of input variables not provided at beginning of testing file ", base.push(P_TESTING_FILE), null);
                                        
                    // Load into an array list each element
                    input = new ArrayList();
                    output = new ArrayList();
                    while(scan.hasNextDouble())
                        {
                        double[] in = new double[numInputs];
                        double out = 0;
                        for(int i = 0; i < numInputs; i++)
                            {
                            if (scan.hasNextDouble())
                                in[i] = scan.nextDouble();
                            else state.output.fatal("Non-normal number of data points in testing file ", base.push(P_TESTING_FILE), null);
                            }
                        if (scan.hasNextDouble())
                            out = scan.nextDouble();
                        else state.output.fatal("Non-normal number of data points in testing file ", base.push(P_TESTING_FILE), null);
                        input.add(in);
                        output.add(new Double(out));
                        }
                                        
                    // dump to arrays
                    len = input.size();
                    testingInputs = new double[len][numInputs];
                    testingOutputs= new double[len];
                    for(int i = 0; i < len; i++)
                        {
                        testingInputs[i] = (double[])(input.get(i));
                        testingOutputs[i] = ((Double)(output.get(i))).doubleValue();
                        }
                    }
                catch (NumberFormatException e)
                    {
                    state.output.fatal("Some tokens in the file were not numbers.");
                    }
                }
            }
        else
            {
            // determine benchmark
            for(int i = 0; i < names.length; i++)
                if (names[i].equals(problem))  // got it
                    { benchmark = i ; break; }
            if (benchmark == -1) // uh oh
                state.output.fatal("Could not find benchmark " + problem, base.push(P_PROBLEM_TYPE), null);
                        
            state.output.message("Doing benchmark " + names[benchmark]);

            try
                {
                trainingInputs = trainPoints(state, benchmark, 0);
                trainingOutputs = new double[trainingInputs.length];
                for(int i = 0 ; i < trainingOutputs.length; i++)
                    trainingOutputs[i] = func(state, trainingInputs[i], benchmark);
                }
            catch (IllegalArgumentException e)
                {
                state.output.fatal("Error in generating training data: " + e.getMessage());
                }

            try
                {
                testingInputs = testPoints(state, benchmark, 0, trainingInputs);
                testingOutputs = new double[testingInputs.length];
                for(int i = 0 ; i < testingOutputs.length; i++)
                    testingOutputs[i] = func(state, testingInputs[i], benchmark);
                }
            catch (IllegalArgumentException e)
                {
                state.output.fatal("Error in generating testing data: " + e.getMessage());
                }

            }
                        
                
        Parameter param = new Parameter("gp.tc.0.fset");  // we assume we have a single tree
        String pval = state.parameters.getString(param, null);
                
        // verify the number of variables match the expected function set
        if (problem == null)  // it's being loaded from file
            {
            boolean found = false;
            String[] vars = fs_vars[trainingInputs[0].length];
            for(int i = 0; i< vars.length; i++)
                if (pval.equals(vars[i])) { found = true; break; }
            if (!found)
                state.output.warning("The number of variables in your problem data (" + trainingInputs[0].length +
                    "does not match the variables found in the function set " + pval + ".  Hope you know what you're doing.",
                    param);
            else state.output.message("Using function set " + pval);
            }
        else
            {
            if (!(pval.equals(fs[benchmark])))  // uh oh
                state.output.warning("The number of variables for the " + names[benchmark] + 
                    " problem (" + trainingInputs[0].length +
                    ") is normally handled by the function set " + fs[benchmark] +
                    " but you are using " + pval + ".  Hope you know what you're doing.  "+ 
                    "To correct this, try adding the parameter gp.tc.0.fset=" + fs[benchmark],
                    param);
            else state.output.message("Using function set " + pval);
            }
        }






    ///// Evaluation.  evaluate(...) uses training cases, and describe(...) uses testing cases


    public void evaluate(EvolutionState state, Individual ind, int subpopulation, int threadnum)
        {
        if (!ind.evaluated)  // don't bother reevaluating
            {
            RegressionData input = (RegressionData)(this.input);

            int hits = 0;
            double sum = 0.0;
            for (int y=0;y<trainingInputs.length;y++)
                {
                currentValue = trainingInputs[y];
                ((GPIndividual)ind).trees[0].child.eval(
                    state,threadnum,input,stack,((GPIndividual)ind),this);

                double error = error(input.x, trainingOutputs[y]);
                                
                // We'll keep the auxillary hits measure for tradition only 
                final double HIT_LEVEL = 0.01;
                if (error <= HIT_LEVEL) hits++; 

                sum += error;              
                }
                
            // the fitness better be KozaFitness!
            KozaFitness f = ((KozaFitness)ind.fitness);
            f.setStandardizedFitness(state, sum);
            f.hits = hits;
            ind.evaluated = true;
            }
        }


    public void describe(EvolutionState state, Individual ind, int subpopulation, int threadnum, int log)
        {
        RegressionData input = (RegressionData)(this.input);

        // we do the testing set here
        
        state.output.println("\n\nPerformance of Best Individual on Testing Set:\n", log);
                
        int hits = 0;
        double sum = 0.0;
        for (int y=0;y<testingInputs.length;y++)
            {
            currentValue = testingInputs[y];
            ((GPIndividual)ind).trees[0].child.eval(
                state,threadnum,input,stack,((GPIndividual)ind),this);

            double error = error(input.x, testingOutputs[y]);
                        
            // We'll keep the auxillary hits measure for tradition only 
            final double HIT_LEVEL = 0.01;
            if (error <= HIT_LEVEL) hits++; 

            sum += error;              
            }
                        
        // the fitness better be KozaFitness!
        KozaFitness f = (KozaFitness)(ind.fitness.clone());     // make a copy, we're just printing it out
        f.setStandardizedFitness(state, sum);
        f.hits = hits;
                
        f.printFitnessForHumans(state, log);
        }
        
    }
